Non-Gaussian Spikes from Chaotic Billiards in Inflation Preheating 
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A new class of non-Gaussian curvature fluctuations Cpr(x) = 5N{xi) arises from the post-inflation 
preheating behaviour of a non-inflaton fleld x»- Its billiard-like chaotic dynamics imprints regular 
log-spaced narrow spikes in the number of preheating e-folds N{xi)- We perform highly accurate 
lattice simulations of SUSY-inspired quartic inflaton and coupling potentials in a separate-universe 
approximation to compute N{xi) as a function of the (nearly homogeneous) initial condition Xi- 
The super-horizon modes of Xi (x) result in positive spiky excursions in (pr and hence negative grav- 
itational potential fluctuations added to the usual sign-independent inflaton-induced perturbations, 
observably manifested in large cosmic structures and as (polarized) temperature CMB cold spots. 

PACS numbers: PACS: 98.80.Cq, CITA-2009-3, astro-ph/0903.3407 



During early universe inflation, vacuum fluctuations 
in all light fields x transformed on super-horizon 
scales into homogeneous and isotropic Gaussian random 
fields. These are fully defined by their power spectra 

^ ^ {\Xk\^)£^, whose magnitude, ~ [i//(27r)]2, is 
related to the Hubble parameter H at expansion fac- 
tor a when Ha first equals the wavenumber k in ques- 
tion. In standard scenarios, metric curvature fluctua- 
tions Sive derived from an inflaton (j) and are nearly 
Gaussian and often nearly scale invariant, with power 
da'^{k) Ri da'^{k) /{2Mpe), in terms of the acceleration 
parameter e = —^^^ and the reduced Planck mass 
Mp = I/^/SttGn [if- Observational constraints on e are 
< 0.03 at the 95% confidence limit from the large angle 
cosmic microwave background radiation (CMB), hence 
is considerably amplified over (jj. Alternative mechanisms 
utilizing light "non-inflatons" x must overcome this e ef- 
fect. Examples are "curvatons" which temporarily dom- 
inate the energy density after inflation [2], and fields 
which induce spatial variations of couplings which modu- 
late the timing of the now-inhomogeneous post-inflation 
(p)reheating [3|], but are not gravitationally dominant. 

The ^ -source first proposed in and studied here is 
nonlinear resonant preheating inducing expansion factor 
variations from the end of inflation at e = 1 to when the 
equation of state settles to w — 1/3. A sample frame- 
work for this has a potential V — jXc/)'^ + \g^<ti^X^i 
with inflation driven by the first and particle creation 
by the second term. If .g^/A is tuned to be ~ C(l), 
the non-infiaton x(x) is light during inflation and accu- 
mulates quantum fluctuations substantially varying on 
scales much greater than the Hubble scale at the end 
of inflation ^ 10"'^Mp. We note that 5^/A=2 cor- 
responds to a SUSY model. On scales shorter than 
H~^, this x(^) is nearly homogeneous defining "sep- 
arate universes" with specified Xi as initial conditions 
for lattice simulations of the fully coupled fields to de- 
termine the nonlinear evolution w, which imprints it- 
self on Ina(xi)- The resulting curvature perturbations 



C,ph (x) are differences in the number of e-folds on uniform 
Hubble (i.e. uniform energy density) time hypersurfaces, 
C, = 5\na\H = SN{xi) 0. If such variations exist, they 
would be in addition to the standard ones. Preheating in 
the model begins with parametric resonance amplifying 
the fiuctuations Xfe(i)e*'^'^ describing vacuum excitations 
of x-particles. The problem can be reduced to a flat 
space-time model by a conformal transformations, Xk(t) 
obeys an oscillator equation with a periodic frequency 
controlled by the background inflaton oscillations 4>{t) 
with (conformal time) period T « 7.416 In the reso- 
nant bands, Xk{t) ~ e'^'"* is unstable, an exponent fik is 
a function of g^/X and k. If g^/X ^ 2, the maximum of 
fik in the first resonant band is located close to A; = 0. 
Hence a nearly homogeneous Xi will be exponentially un- 
stable in each separate universe, and soon becomes en- 
tangled in the complicated mode-mode dynamics of pre- 
heating driven by the back-reaction effects of copiously 
produced x and S(j) particles. To determine SN{x) at 
the part per million level, we need non-linear lattice sim- 
ulations with energy conservation accuracy (EGA) well 
below this. The early attempt by [1] used first-order lat- 
tice simulations with EGA ~ 10~^ to claim an effect at 
the level SN ~ 10~'^, with a dominantly-quadratic form 
SN « /xX? characterized by a constant similar to /nl 
used in non-gaussianity (nG) studies of the CMB sky 0] ■ 
We do not confirm the large quadratic nC of but 
do find a nontrivial SN{xi) with a regular sharp-spiked 
pattern at the observationally interesting ~ 10~^ level, 
as shown in Fig. [1] with a radically different impact on 
the sky than the /nl story. Although our discovery of 
the spiked SN{xi) function was a truly numerical one, a 
posteriori we can explain it by a combination of chaotic 
zero-mode dynamics after te and the abrupt onset of in- 
homogeneous nonlinearity at a preheating time tpr, al- 
lowing us to conjecture for which models it works. 

To accurately compute the very low levels of 6N{xi), 
especially since we are in finding and validation mode for 
such tiny effects, we needed to go well beyond the EGA 
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FIG. 1: The structure of SN{xi) on uniform Hubble hypersurfaces probed with ~ 10* lattice simulations from the end-of- 
inflation through the end-of-preheating for varying homogeneous Xi initial conditions, for /\ — 2. The periods n/ioT in 
\nxi are marked by the large green circles, and the higher harmonics (revealed by the Fourier analysis) by smaller green 
circles. These locate the spikes in 5N{xi)- The effective field (-FNLlXft + X>h) marginalized over high spatial frequencies with 
cthf=7 X 10~^Mp (vertical line) yields the solid curve. A quadratic fit, fxiXb + X>h)'^, is also shown. An issue for our Hubble 
patch is whether the ultra-large scale x>h is large enough that the large scale structure fiuctuations about it, ±ab<h, encompass 
smoothed peaks in field space, or not. A typical value for a6<h is ~ 3 x 10~^Mp. 



practically achievable in second-order preheating codes 
(typically ~ lO'^ for LATTICEEASY ^ and - 10"^ for 
DEFROST [9]). We developed a new high-order symplec- 
tic PDE solver with adaptive time steps, which can reach 
machine precision levels (EGA - IQ-^^!) to address this 
problem. In retrospect, we find that EGA of 10~^ would 
suffice, which is achievable with shorter timesteps in DE- 
FROST. The lattice simulations begin at expansion fac- 
tor fle when e = 1; we have shown that variations of the 
start do not effect results. The calculations are stopped 
at a uniform Hubble value Hf when the w ~ p/ pis nearly 
1/3, the radiation-dominated value. Although the aver- 
age w is ~ 1/3, it has small oscillations during preheat- 
ing, leading to a fluctuations to average over to determine 
N{xi) for use in the curvature formula. To deal with this, 
we used an accurate extrapolation template for the aver- 
aged a{t) or, for superb accuracy, a Kaiser window filter 
to suppress high frequency oscillations in a{t) by a factor 
10~^. With either, we showed that provided we calcu- 
late for 5-6 e-folds after Oe, the specific Hf is not impor- 
tant. An accuracy test was to show 6N is effectively zero 
(<C 10^^) and not modulated by Xi for 9^/^ out of the 
resonant band, e.g., at g^ /X = 1 and 3, for which /io = 0. 
With the symplectic code, we ran a 256'^ lattice with a 
box size 64/ He to check that our conclusions derived us- 
ing a very large number of lower resolution simulations 
using DEFROST to build up statistics in Xi are accu- 
rate. The essential effects can indeed be captured with 
lower resolution and box size. Fig. [1] showing a spiky 
pattern in 6N{xi) for g^/X = 2 was produced with 11563 
32'^-simulations with box size 20/ H^. The amplitude of 
spikes increases as Xi increases. The spikiest pattern is 
at g^/X = 1.875, with the spikes broadening away from 
that, finally disappearing at the g^/X = 1 and 3 borders. 



For Fig. [1] we explored a large range in input Xi/Mp, 
from 10"^ up to 2 X 10"^. Xi(x) = X<h + X>h has 
a sub-horizon contribution x<h from eigenmodes with 
wavenumbers between iJeie to the current hori- 

zon scale, kh ~ Hhau, and a super-horizon contribution 
X>h with waves from ku to a fcmm whose value will de- 
pend upon the inflation model. The corresponding vari- 



ances are cr; 



fe, '^^x and ^> 
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> 100 e-folds, da^ is nearly constant for 1 < g^/X < 2, 
drops substantially as a — > Oe for g^ /X — cj)'^ as ex- 
pected), and actually rises for g^ /X = 3. The least num- 
ber of e-folds Ina/i/fle must exceed 50 — 60, and since 
H > He during inflation, cr|^ > Inke/kh [i?e/(27r)]^ 
gives a x<h enhanced over He by ~ 1/55. The super- 
horizon power, cr?.^ > Inkh/kmin [He/i^n)]'^ is also log- 
enhanced, and considerably so in our illustrative Xcj)'^ ex- 
ample. Thus, the log factors give larger Xii including a 
X>h random number which is nearly constant within our 
Hubble patch, but has a ~ ±cr>h patch-to-patch "cosmic 
variance" . 

We now show how the features of SN(xi) can be 
understood qualitatively from trajectories in the two- 
dimensional space of homogeneous modes {(j){t) , x{t)) ■ 
The excited inhomogeneous degrees of freedom do back- 
react on these /c « modes, but only later in the evo- 
lution, e.g., at t ~ lOT for g'^/X = 2. The {4>{t) , x{t)) 
space is effectively bounded by potential energy barriers 
jX(j)'^ + ^g'^(l)'^X^ = const, as shown in Fig.O Initially the 
trajectories oscillate mostly in the (fi direction, with only 
very small initial amplitudes in Xj a-s illustrated in the 
insets in Fig. [2l These oscillations are akin to billiard 
motions between the potential walls. Precession of the 
initial oscillations causes the x amplitude to grow expo- 
nentially in a chaotic manner: x(i) = Xi e^*, where A is 
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FIG. 2: Billiard trajectories of the A; = modes (j){t) and 
within the A(/!)*/4 + g'^x^'t'^/^ potential well. Upper panel 
has a "no-spike" initial value X' = 3.6 x W~^Mp, and the 
bottom panel has a "spike" X' = 3.9 x 10~^ Mp. The solid 
curves are the (fuzzy) potential wall without the inclusion of 
mass terms induced by field nonlinearities; the pale green and 
brighter green border areas include the induced masses at the 
instances t = 10. ST and asymptotically. Thin and thick parts 
of the trajectories denote before after t — 9.7T, and up to 
11. Sr, with the circle on it at t = 10. ST. The inserts in the 
left upper corners of the panels show the first several periods 
of linear oscillations. 

the Lyapunov exponent. This gives us new insight on the 
parametric resonant k = solution in terms of the Lya- 
punov instability. For g^/X = 2, we find A = 0.235. 
This conjecture also works for the cases in which fc = 
is not in the resonant band; e.^.jfor g^/X = 3 we have 
fiQ — and do indeed find periodic trajectories so A = 0. 

As X grows beyond the linear regime, the bouncing bil- 
liard experiences the negative curvature of the potential 
walls and a bifurcation of the trajectories occurs, with a 
few entering the arms in between walls, and most do not. 
If there were only homogeneous modes the impact of this 
would be temporary because eventually all trajectories 
would be chaotically mixed. However, the excitation of 
the inhomogeneous modes result in exponential growth 
of (x^) and ((50^), as e^'^** and e'''^*' respectively, where 
/i* is an effective resonant exponent [5]. These induce 
enhanced effective masses in the fields, abruptly chang- 
ing the potential, with the arms in between walls closing 
exponentially quickly, as shown in Fig. [2l The trajec- 
tories which happen to get into the arms before arm- 
closure evolve very different from those which never get 



into the arms. The billiard picture breaks down when 
the gradient terms in x and S(j) occur, at arm-closure 
time, but the bifurcation determined by the initial con- 
ditions at the linear stage has already happened. The 
two pre-closure classes of in-arm and not-in-arm trajec- 
tories result in different equations of state, and hence a 
Xi-dependent SN{xi)- In-arm trajectories experience ki- 
netic energy kicks from the closing arms, which translates 
into a transient decrease in w, inducing a jump in SN. 
The no-spike trajectories of the upper panel of Fig. [2] are 
much more numerous than the rarer sort in the lower 
panel corresponds that give spikes. The billiard picture 
predicts the spike pattern as one periodic in Inx which 
works extremely well: the same "spiky" trajectory la- 
beled by Xi is repeated for initial values Xi e"'^""^ for 
integer n. The origin of the higher harmonics is more 
complex. Our spike pattern formula works very well for 
other values of g^/X, with the requisite no{g'^/X). 

Let us denote Cpr by ^NL(Xi(x))- These are added 
to the conventional inflaton curvature fluctuations Crt>- 
Ctot(x) « C0(x) + Fnl(x(x)) - (-Fnl). (The ensemble 
average (Fnl) over all values of x, is subtracted so (tot 
fluctuates about zero.) Structure in Xi and hence in Fnl, 
exists on a vast range of resolution scales, from kmin 
through kh to fee. Observed large scale structure (LSS), 
as probed by the CMB, redshift surveys, weak lensing, 
and rare event abundances such as of clusters probe ~10 
e-folds in In A; below k^^, and about the same number of 
e-folds just below this can be probed with the more un- 
certain astrophysical observables involving earlier stages 
in the nonlinear hierarchy, galaxies, dwarves, and the 
"flrst stars" within them. The superhorizon scales be- 
yond fc^^ have no impact on ^^-structure, but do have 
a large impact on _Fnl through the specific value x>h 
built from waves fcmm < k < k^. To explore the astro- 
physical consequences of i^NL, we marginalize over high 
frequency components xhf of x = Xhf + Xlf to form the 
conditional non-Gaussian "effective field", (-FnlIxlf), 
with ~40-50 e-folds of "short distance" substructure 
in i^NL filtered out. The mean-field, (FnlIxlf) — 
J -Fnl(x)-P(xIxlf)o?x, is a Gaussian-smoothing of -Fnl 
in field-space, via P(xIxlf) - '^^PKx^XLF)V(2cTgip)] ^ ^^^^ 

variance cr^p = (xhf)- Using it will give reliable LSS in- 
ferences if the Fourier transform of the deviation Ai^NL = 
Fnl — (-FnlIxlf) is smafl for k~^ 3> fc^p, the xlf filter 
scale. Since xlf = X6 + X>h contains a spatially varying 
part Xb built from waves with kh < k < /clf and a con- 
stant superhorizon part x>h- Which aspects of the spiky 
patterns of Fig. [T] that would be realized in our Hub- 
ble patch is quite dependent on the luck of our x>/i-draw 
from a Gaussian distribution with variance ct?,^ = iXyh)^ 
in particular whether x>h is near a smoothed-peak, or 
small ^ cthf- 

The analytic model, i^NL = X]p *5xp [— (x — 
Xp)^/(27p)], approximates each spike of Fig. [1] with 
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FIG. 3: Realizations of the nG map (-FkLlXi)) on the CMB 
sky. Top left shows a scale-invariant Gaussian random field 
realization Xbid-ydecq) in direction q on a sphere at the co- 
moving distance to photon decoupling, d-ydec Top right shows 
the action of (i^NLlXft) on it, using our Gaussian-line-profile 
approximation with 2 peaks at Xp ~ ^VpCJb, for Vp = 3.5. 
Middle left shows the map convolved with a CMB transfer 
function, and smoothed on a 1° scale, right with Up = 4.5; 
both show "cold spot" intermittency. 

a Gaussian "line profile" of width 7p centered on Xpi 
with peak amplitude Fp and integrated line strength 
Up = Fp^/2TT^p. The X — > — X symmetry means that 
for each peak at Xp there is a mirror peak at — Xp of 
the same strength. The conditional n-point moments 
of i^NL, (nr=i ^NL(xi)lxLF) are then easily computable 
Gaussians, with quadratic terms in Xi"Xp linked through 
the HF 2-point function of xhf- The mean field has 
n = 1: 

p '^PHF 

Here cr^jjp=7p -I- cr^p is « cr^p for the typical 7p we find 
and [/p=^exp{-^}, with Xp>h--Xp-X>h- The 
Xhf correlation function dependence of the fluctuation 
variance (A_F'nl(xi)AFnl(x2)|xlf) precludes an ana- 
lytic harmonic analysis, but we have investigated this 
numerically: e.g., for a 1-spike version, this i^NL variance 
is 1% of (FnlIxlf)^ if it is smoothed on ~ fc^p, and 
^ 0.1% with ~ lOfc^p smoothing. (We have checked 
the effective field works even better for Fnl quadratic 
and exponential in x-) The HF structure is very rel- 
evant to first-object formation, but what actually hap- 
pens will be model-dependent. However, HF will not 
impinge upon LSS observables, since they are convolved 
with experimentally-determined or theoretically-imposed 
windows ^ in scale. 

A fundamental character of the resonant mechanism 
is the delay of in-arm preheating completion, translat- 
ing into positive large excursions in 5\na\H. The asso- 
ciated perturbed Newtonian gravitational potential $n 
is negative. The CMB sky temperature T{q), an an- 
gular function of the CMB photon direction q towards 
us, is a projected image of various sources whose 3D 
Fourier transforms involve various form factors J-"(|k|) 
times <l>N(k, to)/3- The dominant T terms are from two 
CMB decoupling effects, and one late-time effect: the 
combined "naive Sachs- Wolfe (NSW) effect plus photon 



compression-rarefaction"; the Doppler effect from flow- 
ing electrons; and the integrated Sachs- Wolfe (ISW) ef- 
fect with a fc-dependent time-integral of 6<i>N/$N 
The upper panels in Fig. [3] correspond to J- = 1 and the 
lower panels are convolved with a CMB transfer func- 
tion and smoothing on 1°, appropriate if the CMB sky is 
a direct map of the photon decoupling surface, ignoring 
its fuzziness (valid for the 1° smoothing) and a correct 
implementation of the ISW effect, which will affect the 
largest scales. However, the essential intermittent cold 
spot nature following from negative will persist. Such 
cold spots will be polarized just as those deriving from 
are, and have a CMB spectrum, a prediction for the ~ 6° 
COBE/WMAP cold spot 0] if it is a smoothed remnant 
of resonant preheating. 

Another non-intermittent nG regime is that of ex- 
treme "line blending": Fig. [T] shows IxXIf/^'^p with 
f^= 2 X 10^ is a good fit to (FnlIxlf) below ~ cthf- 
(Such a quadratic form also follows from Eq. For 
X>h ;S CTHF, we will get a power law in Xfe, PxXb/Mp + 
fxXb/Mp: with f3x='^fxX>h/Mp. The conventional 
WMAP -9 < /nl < 111 constraints use - C</.+/nlC^- 
In our case, /xX^ is uncorrected with so the constraint 
on /x will be considerably relaxed, as long as the linear 
PxXb term is subdominant to as is expected. 

Further exploration is needed of how spiked Fnl from 
resonant preheating may arise in more general inflation 
models and on the rich nG impact and observable con- 
straints of such Fnl on short and long cosmic scales. For 
now, we note SUSY models provide light non-inflaton 
fields of the sort we need, and future CMB experiments 
could test whether nG cold spots are polarized. We 
thank A. Chambers, G. Felder, E. Komatsu, A. Linde 
and A. Rajantie for discussions. 
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